AFRL-RV-PS- 

TR-2016-0161 


AFRL-RV-PS- 

TR-2016-0161 


ATOMISTIC-  AND  MESO-SCALE 
COMPUTATIONAL  SIMULATIONS  FOR 
DEVELOPING  MULTI-TIMESCALE  THEORY  FOR 
RADIATION  DEGRADATION  IN  ELECTRONIC  AND 
OPTOELECTRONIC  DEVICES 


Fei  Gao 


University  of  Michigan 
2355  Bonisteel  Blvd. 

Ann  Arbor,  Michigan  48109-2014 


13  Feb  2017 


Final  Report 


APPROVED  FOR  PUBLIC  RELEASE;  DISTRIBUTION  IS  UNLIMITED. 


AIR  FORCE  RESEARCH  LABORATORY 

Space  Vehicles  Directorate 

3550  Aberdeen  Ave  SE 

AIR  FORCE  MATERIEL  COMMAND 

KIRTLAND  AIR  FORCE  BASE,  NM  87117-5776 


DTIC  COPY 


NOTICE  AND  SIGNATURE  PAGE 

Using  Government  drawings,  specifications,  or  other  data  included  in  this  document  for 
any  purpose  other  than  Government  procurement  does  not  in  any  way  obligate  the  U.S. 
Government.  The  fact  that  the  Government  formulated  or  supplied  the  drawings, 
specifications,  or  other  data  does  not  license  the  holder  or  any  other  person  or  corporation; 
or  convey  any  rights  or  permission  to  manufacture,  use,  or  sell  any  patented  invention  that 
may  relate  to  them. 


This  report  is  the  result  of  contracted  fundamental  research  deemed  exempt  from  public  affairs 
security  and  policy  review  in  accordance  with  SAF/AQR  memorandum  dated  10  Dec  08  and 
AFRL/CA  policy  clarification  memorandum  dated  16  Jan  09.  This  report  is  available  to  the  general 
public,  including  foreign  nationals.  Copies  may  be  obtained  from  the  Defense  Technical 
Information  Center  (DTIC)  (http://www.dtic.mil). 


AFRL-RV -PS -TR-20 16-0161  HAS  BEEN  REVIEWED  AND  IS  APPROVED  FOR 
PUBLICATION  IN  ACCORDANCE  WITH  ASSIGNED  DISTRIBUTION  STATEMENT. 


//SIGNED// 

DAVID  CARDIMONA 

Technical  Advisor,  Space  Based  Advanced  Sensing 
and  Protection 


//SIGNED// 

JOHN  BEAUCHEMIN 

Chief  Engineer,  Spacecraft  Technology  Division 
Space  Vehicles  Directorate 


//SIGNED// 

DAVID  CARDIMONA 
Program  Manager 


This  report  is  published  in  the  interest  of  scientific  and  technical  information  exchange,  and  its 
publication  does  not  constitute  the  Government’s  approval  or  disapproval  of  its  ideas  or  findings. 


Approved  for  public  release;  distribution  is  unlimited. 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMB  No.  0704-0188 


Public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and  maintaining  the 
data  needed,  and  completing  and  reviewing  this  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information,  including  suggestions  for  reducing 
this  burden  to  Department  of  Defense,  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports  (0704-0188),  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington,  VA  22202- 
4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  any  penalty  for  failing  to  comply  with  a  collection  of  information  if  it  does  not  display  a  currently 
valid  OMB  control  number.  PLEASE  DO  NOT  RETURN  YOUR  FORM  TO  THE  ABOVE  ADDRESS. 


1 .  REPORT  DATE  (DD-MM-YY)  2.  REPORT  TYPE  3.  DATES  COVERED  (From  -  To) 

13-02-2017  Final  Report  18  Aug  2015  -  7  Dec  2016 


4.  TITLE  AND  SUBTITLE  5a.  CONTRACT  NUMBER 

Atomistic-  and  Meso-Scale  Computational  Simulations  for  Developing  Multi-Timescale 

Theory  for  Radiation  Degradation  in  Electronic  and  Optoelectronic  Devices  FA9453-15-1-0084 

5b.  GRANT  NUMBER 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

University  of  Michigan 
2355  Bonisteel  Blvd. 

Ann  Arbor,  Michigan  48109-2014 


5c.  PROGRAM  ELEMENT  NUMBER 

6260 IF 


5d.  PROJECT  NUMBER 

4846 


5e.  TASK  NUMBER 

PPM00015486 


5f.  WORK  UNIT  NUMBER 

EF126376 


8.  PERFORMING  ORGANIZATION  REPORT 
NUMBER 


9.  SPONSORING  /  MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 

Air  Force  Research  Laboratory 
Space  Vehicles  Directorate 
3550  Aberdeen  Ave.,  SE 
Kirtland  AFB,  NM  87117-5776 


12.  DISTRIBUTION  /  AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  is  unlimited. 


10.  SPONSOR/MONITOR’S  ACRONYM(S) 

AFRL/RVSW 


11.  SPONSOR/MONITOR’S  REPORT 
NUMBER(S) 

AFRL-R  V-PS  -TR-20 1 6-0 1 6 1 


14.  ABSTRACT 

This  project  will  support  a  research  effort  at  the  University  of  Michigan  to  simulate  atomistic-  and  meso-scale  behavior  of  defect  evolutions 
in  compound  semiconductors,  including  ultrafast  displacement  cascade,  intermediate  defect  stabilization  and  cluster  formation,  as  well  as 
slow  defect  reaction  and  migration.  The  fundamental  mechanisms  and  knowledge  gained  from  atomic-  and  meso-scale  simulations  will  be 
input  into  rate-diffusion  theory  as  initial  conditions  to  calculate  the  steady-state  distribution  of  point  defects  in  a  mesoscopic  layered- 
structured  system,  thus  allowing  the  development  of  a  multi-timescale  theory  to  study  radiation  degradation  in  electronic  and  optoelectronic 
devices. 


15.  SUBJECT  TERMS 

barrier,  detector,  infrared,  unipolar 


16.  SECURITY  CLASSIFICATION  OF: 

a.  REPORT 

Unclassified 

b.  ABSTRACT 

Unclassified 

c.  THIS  PAGE 

Unclassified 

17.  LIMITATION 

18.  NUMBER 

OF  ABSTRACT 

OF  PAGES 

Unlimited 

30 

19a.  NAME  OF  RESPONSIBLE  PERSON 

David  Cardimona 


19b.  TELEPHONE  NUMBER  (include  area 
code) 


Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std.  239.18 


(This  page  intentionally  left  blank) 


Approved  for  public  release;  distribution  is  unlimited. 


Table  of  Contents 


1.  SUMMARY - - - - — - - - - - . — - .  1 

2.  INTRODUCTION - - - - - - - - - 1 

3.  METHODS,  ASSUMPTIONS,  AND  PROCEDURES - - - - - 2 

3.1.  Molecular  Dynamics  Model  of  Displacement  Energy - 2 

3.2.  Molecular  Dynamics  Model  of  Defect  Generation  — - - - —  4 

3.3.  Interatomic  Potentials . . — . .  4 

3.4.  NIEL  Model -  5 

4.  RESULTS  AND  DISCUSSIONS - - — -  6 

4.1.  Displacement  Energies - 6 

4.1.1.  Unit  Stereographic  Triangles - 7 

4.1.2.  Threshold  Displacement  Energy  Contour  Maps - - -  7 

4.1.3.  Overview  of  Ed - - — - - - -  9 

4.2.  Displacement  Cascades  and  Effective  NIEL - - - 9 

4.2.1.  Defect  Generation - - 9 

4.2.2.  Effective  NIEL - - - - -  14 

5.  CONCLUSIONS - - - - - - - - -  16 

REFERENCES - - - - - - - 18 

LIST  OF  ACRONYMS - - - - - . .  21 


Approved  for  public  release;  distribution  is  unlimited. 


LIST  OF  FIGURES 


Figure  Page 

1.  (a)  GaAs  crystal  structure;  (b)  illustration  of  stereographic  triangles  on  a  reference  sphere.... 3 

2.  (a)  The  displacement  threshold  energy  around  the  sides  of  the  unit  stereographic  triangles 
for  Ga;  (b)  the  displacement  threshold  energy  around  the  side  of  the  unit  stereographic  triangles 


for  As . 7 

3.  (a)  Threshold  energy  map  for  Ga;  (b)  threshold  energy  map  for  As . 8 


4.  Atomistic  configurations  of  a  typical  10  keV  cascade  in  GaAs  at  (a)  the  peak  damage  and 
(b)  the  final  damage  state,  where  only  defects  are  shown.  The  defects  are  identified  by  size  and 


color,  as  indicated  by  different  spheres . 10 

5.  The  number  of  displaced  atoms,  including  Ga  and  As  displacements  and  antisites,  as  a 
function  of  time  for  the  10  keV  cascade  shown  in  Fig.  4 . 11 


6.  The  final  defect  states  of  (a)  2  keV,  (b)  5  keV  and  (c)  20  keV  cascade  in  GaAs  from  which  can 


be  clearly  seen  the  transition  from  a  single  pocket  of  atomic  displacements  to  multiple 
subcascades . 12 

7.  The  number  of  surviving  defects  produced  in  the  cascades  as  a  function  of  PKA  energy, 

where  Ga  and  As  interstitials,  as  well  as  antisite  defeats,  are  presented  separately . 13 

8.  The  total  defects  of  the  final  damages  as  a  function  of  the  energy  density  deposited  by  a  PKA, 

where  the  solid  line  represents  the  fits  to  MD  data . 14 


9.  The  NIEL  obtained  based  on  the  MD  results  for  several  different  incident  particles  in  GaAs, 
where  the  values  of  the  NIEL  calculated  by  an  analytical  approach  are  imposed  for 
comparison . 16 


Approved  for  public  release;  distribution  is  unlimited. 


ACKNOWLEDGMENTS 


This  material  is  based  on  research  sponsored  by  Air  Force  Research  Laboratory  under 
agreement  number  FA9453-15-1-0084.  The  U.S.  Government  is  authorized  to  reproduce  and 
distribute  reprints  for  Governmental  purposes  notwithstanding  any  copyright  notation  thereon. 


DISCLAIMER 

The  views  and  conclusions  contained  herein  are  those  of  the  authors  and  should  not  be 
interpreted  as  necessarily  representing  the  official  policies  or  endorsements,  either  expressed  or 
implied,  of  Air  Force  Research  Laboratory  or  the  U.S.  Government. 


Approved  for  public  release;  distribution  is  unlimited. 


(This  page  intentionally  left  blank) 


iv 


Approved  for  public  release;  distribution  is  unlimited. 


1.  SUMMARY 


This  project  will  support  a  research  effort  at  the  University  of  Michigan  to  simulate 
atomistic-  and  meso-scale  behavior  of  defect  evolutions  in  compound  semiconductors,  including 
ultrafast  displacement  cascade,  intennediate  defect  stabilization  and  cluster  fonnation,  as  well  as 
slow  defect  reaction  and  migration.  The  fundamental  mechanisms  and  knowledge  gained  from 
atomic-  and  meso-scale  simulations  will  be  input  into  rate-diffusion  theory  as  initial  conditions  to 
calculate  the  steady-state  distribution  of  point  defects  in  a  mesoscopic  layered-structured  system, 
thus  allowing  the  development  of  a  multi-timescale  theory  to  study  radiation  degradation  in 
electronic  and  optoelectronic  devices.  The  long-term  goal  is  developing  a  fundamental 
understanding  of  defects  and  defect  processes  in  compound  semiconductors,  including 
defect/property  relationships,  the  effects  of  defects  on  transport  processes,  the  aggregation  of 
defects  to  fonn  complex  nanostructures,  and  the  development  of  predictive  models  of  behavior. 


2.  INTRODUCTION 

It  is  well  known  that  in  a  perfect  crystal,  the  continuous  free-electron  states  are  quantized 
into  many  Bloch  bands  separated  by  energy  gaps.  These  Bloch  electrons  move  freely  inside  the 
crystal  with  an  effective  mass  different  from  the  free-electron  mass  [1].  In  the  presence  of  defects, 
however,  the  field-driven  current  flow  of  Bloch  electrons  in  the  perfect  crystal  will  be  scattered 
locally  by  these  defects,  leading  to  a  reduced  electron  mobility.  Also,  the  photo-excited  electron 
lifetime,  due  to  non-radiative  recombination  with  defects,  has  been  proven  to  be  a  key  factor 
affecting  the  sensitivity  or  the  perfonnance  of  optoelectronic  devices  (e.g.,  photo-detectors  and 
light-emitting  diodes)  [2],  The  dangling  bonds  attached  to  the  point  defects  may  capture  extra 
electrons  to  form  charged  defects.  In  this  case,  the  positively-charged  holes  in  the  system  will  be 
trapped  to  produce  a  strong  space-charge  field,  while  the  negatively-charged  electrons  may 
generate  the  so-called  1/f-current  noise  in  their  bumpy  motions  due  to  the  presence  of  many 
potential  minima  and  maxima  from  randomly-distributed  charged  defects.  Therefore, 
understanding  of  defect  production,  stabilization,  clustering,  migration  and  interaction  with 
microstructural  features  is  crucial  for  developing  a  multi-scale  theory  to  explore  radiation 
degradation  in  electronic  and  optoelectronic  devices.  At  the  most  fundamental  level,  molecular 
dynamics  can  be  used  to  study  defect  production,  migration  and  interaction,  while  kinetic  Monte 
Carlo  methods  can  be  employed  to  simulate  defect  evolution  and  their  spatial  distribution.  The 
researchers  at  University  of  Michigan  have  extensive  experience  performing  atomic-  and  meso- 
scale  simulations  of  defect  production  and  evolution  in  semiconductors. 

The  current  project  focuses  on  one  particular  semiconductor,  namely  gallium  arsenide 
(GaAs),  because  it  has  received  considerable  attention  due  to  its  potential  electronic  applications, 
such  as  GaAs-based  metal  semiconductor  field  effect  transistors  and  logic  gates,  near-infrared 
imaging  devices  and  photovoltaic  nodules  [3],  As  a  means  for  fabricating  compound 
semiconductors,  the  interest  in  using  ion  implantation  continuously  increases,  which  inevitably 
introduces  defects.  The  diffusion  and  accumulation  of  these  defects  are  critical  factors  controlling 
the  dynamics  of  ion  implantation  processes  and  device  degradation,  as  well  as  affecting  charge 
compensation,  minority-carrier  lifetimes  and  luminescence  efficiencies.  Recently,  use  of  GaAs  in 
high-power  space-energy  systems  and  special  space-probe  applications  has  been  proposed. 
However,  space  radiation  damage  to  the  GaAs  may  be  a  limiting  factor  on  interplanetary  missions 
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unless  sufficient  shielding  is  provided  to  keep  damage  levels  under  acceptable  limits. 
Consequently,  radiation  damage  studies  have  been  made  experimentally  on  the  effects  of  electron 
[4],  and  proton  and  neutron  irradiation  [5],  including  defect  production  and  annealing,  as  well  as 
effects  on  the  perfonnance  of  GaAs  devices.  On  the  other  hand,  understanding  the  basics  of  ion- 
solid  interaction  and  irradiation  damage  has  led  to  significant  developments  in  state-of-the-art 
atomic-level,  kinetic  Monte  Carlo  and  meso-scale  simulations,  and  these  simulations  have 
dramatically  advanced  the  knowledge  of  defects  and  defect  processes  in  a  number  of  materials, 
ranging  from  metals  [6,7]  to  semiconductors  [8]  to  ceramics  [9,10].  Recently,  large-scale  ab  initio 
and  classical  molecular  dynamics  (MD)  methods  have  been  developed  for  the  study  of  radiation 
damage  in  semiconductors  [11,12],  and  these  methods  are  used  to  explore  the  number  of 
displacements  produced,  defect  clustering  and  disordering,  as  well  as  the  effects  of  charge  transfer 
and  charge-density  redistribution  on  the  dynamics  and  ultimate  charge-state  of  defect  fonnation. 
These  simulations  have  demonstrated  some  nonlinear  effects  taking  place  at  low  and  high  energies, 
which  can  greatly  modify  the  number  of  displacements  as  predicted  by  the  simplified  Kinchin- 
Pease  model. 

It  has  been  long  realized  that  the  observed  radiation  damage  to  electronic  components  made 
from  semiconductors  or  semiconducting  compounds  is  proportional  to  the  non-ionizing  energy 
loss  (NIEL)  [13],  Consequently,  the  experimental  use  of  NIEL  for  correlating  proton-induced 
displacement  damage  in  semiconductor  devices  has  been  widely  applied  over  the  past  decade, 
which  has  proven  useful  in  the  study  of  both  Si  and  GaAs  [14,15].  Moreover,  radiation  damage  in 
other  semiconductors  such  as  diamond  and  SiC  [16],  and  many  III-V  compounds,  such  as  InP  [17], 
has  been  investigated  using  the  NIEL  concept.  Previously,  several  models  have  been  developed  to 
estimate  the  NIEL  in  GaAs  using  a  Monte  Carlo  charged  particle  transport  code  [18]  or  based  on 
an  empirically  determined  damage  efficiency  [19].  Although  these  models  provide  significant 
insights  into  the  correlation  between  the  NIEL  and  displacement  damage,  there  remain 
discrepancies  between  proton  NIEL  calculations  and  experimental  measurements  in  GaAs  devices. 
Determining  accurate  damage  levels  and  improving  the  NIEL  model  have  been  an  ongoing  focus 
of  the  Space  Radiation  Effects  Community  (SREC).  In  light  of  this,  we  have  employed  a  molecular 
dynamics  method  to  simulate  defect  production,  defect  clustering  and  disordering,  as  well  as  initial 
kinetic  energy  density  in  GaAs.  In  these  simulations,  the  damage  state  following  each  cascade  was 
characterized  and  analyzed  for  defect  distribution  and  amorphization.  Based  on  MD  results,  we 
have  developed  a  new  energy  partition  function  model,  which  was  then  used  to  determine  the  NIEL 
in  GaAs.  The  present  approach  provides  a  general  method  for  accurately  calculating  the  NIEL  for 
fast  neutrons,  protons  and  heavier  particles. 


3.  METHODS,  ASSUMPTIONS,  AND  PROCEDURES 
3.1.  Molecular  Dynamics  Model  of  Displacement  Energy 

A  modified  version  of  the  Molecular  Dynamics  (MD)  code-MOLDY  [20]  was  employed  to 
determine  the  threshold  displacement  energies  along  different  recoil  directions.  A  simulation  box 
with  a  size  of  10x10x10  unit  cells  (consisting  of  8,000  atoms)  was  created,  with  a  zinc-blende 
structure  for  GaAs,  as  shown  in  Fig.  1(a).  Periodic  boundary  conditions  were  applied  to  the  system. 
During  the  simulation,  the  computational  block  was  held  at  a  temperature  at  0  K,  and  equilibrated 
for  5  ps  prior  to  a  recoil  event.  To  simulate  a  recoil  event,  a  recoiling  atom  was  created  by  defining 
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either  Ga  or  As  as  a  primary  knock-on  atom  (PKA),  which  was  chosen  close  to  the  center  of 
simulation  block  to  prevent  the  PKA  from  crossing  the  boundaries.  Furthermore,  the  PKA  was 
given  a  kinetic  energy  (starting  at  6  eV)  along  a  specific  crystallographic  direction  and  the  dynamic 
process  was  evolved  for  5000  MD  time  steps  (about  5ps  in  real  time).  If  the  amount  of  energy 
given  was  enough  to  overcome  the  threshold  displacement  energy,  the  PKA  was  initiated  by 
displacing  it  from  its  original  site  to  generate  a  stable  Frenkel  pair  or  anion  anti-site  at  the  end  of 
the  simulation.  On  the  other  hand,  if  a  stable  defect  configuration  was  not  observed,  the  initial 
kinetic  energy  was  increased  by  2  eV,  and  so  Ed,  quoted  here  is  ±1  eV. 


Fig.  1  (a)  GaAs  crystal  structure;  (b)  illustration  of  stereographic  triangles  on  a  reference  sphere. 


Consideration  of  crystallographic  geometry  has  been  proven  to  be  a  distinctive  factor  in 
exploring  the  displacement  damage  [21].  In  this  study,  we  use  the  Mercator  method  [22]  to 
represent  the  stereographic  projection  of  the  crystallographic  direction  of  the  zinc-blende  crystal, 
which  can  be  thought  to  be  equivalent  to  the  diamond  structure.  Each  point  on  the  reference  sphere 
in  Fig.  1(b)  stands  for  a  direction  projection,  where  the  sphere  center  coincides  with  the  origin  of 
the  crystal  cube.  For  zinc-blende,  a  standard  unit  stereographic  triangle  on  the  sphere  is 
representative  of  the  area  in  which  all  possible  recoil  directions  can  take  place,  with  the  remaining 
directions  known  to  exhibit  the  same  behavior  [22].  However,  for  a  compound  semiconductor,  e.g. 
GaAs,  which  has  the  nearest  neighbor  of  unlike  atoms,  two  triangles  are  required  (with  the  side  of 
[100]-[1 10]-[1 1 1]  and  [100]-  [110]  -[Ill],  respectively)  to  investigate  the  displacement  threshold 
energy  for  each  atom  type.  For  example,  an  As  atom  moving  along  the  [111]  direction,  as  indicated 
by  an  orange  arrow  in  Fig.  1(a),  encounters  a  Ga  atom  at  the  distance  of  7.34  A,  but  at  a  distance 
of  2.45  A  when  it  moves  along  the  [111]  direction.  This  should  make  the  displacement  of  As  atoms 
difficult  along  [111]  and  much  easier  along  [111].  This  fact  has  been  illustrated  in  Fig.  1(b),  with 
red  shaded  part  corresponding  to  the  triangle  in  positive  direction  while  the  blue  one  is  for  negative. 
Thus,  more  than  3600  directions  in  total  were  simulated  for  both  species  as  a  PKA  indicated  by 
the  angles  in  Mercator  map. 
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3.2.  Molecular  Dynamics  Model  of  Defect  Generation 

To  simulate  displacement  cascades,  periodic  boundary  conditions  were  also  applied  in  three 
directions,  and  the  lattice  temperature  was  initially  selected  to  be  100  K,  which  mimics  a  thermal 
energy  just  about  the  zero-point  value.  All  the  simulations  are  carried  out  with  constant  volume 
and  constant  number  of  atoms.  The  MD-cell  size  is  selected  to  avoid  overlap  of  a  displacement 
cascade  with  itself  due  to  the  periodicity.  For  a  displacement  cascade  of  20  keV,  a  crystal  of 
1,000,000  atoms  (50x50x50  unit  cells)  is  employed  for  simulation.  Before  the  primary  recoil  is 
initiated,  the  perfect  crystal  block  is  allowed  to  equilibrate  for  10  ps  at  100  K  to  achieve  an 
equilibrium  phonon  state.  The  cascades  are  initiated  by  giving  a  primary  knock-on  atom  (PKA)  a 
kinetic  energy  along  a  randomly-selected  direction.  In  the  present  study,  only  Ga  PKA  events 
are  simulated,  but  the  simulations  can  be  extended  to  As  PKAs.  Furthermore,  the  collision  cascade 
produces  secondary  As  knock-on  atoms,  which  contribute  to  the  overall  damage  production.  The 
evolution  of  each  cascade  is  examined  for  about  20  ps  following  the  primary  knock-on  event,  after 
which  the  number  of  defects  shows  no  significant  change.  We  have  simulated  the  displacement 
cascades  in  the  energy  range  from  500  eV  to  20  keV,  and  20  PKA  events  are  carried  out  at  each 
energy  to  obtain  meaningful  statistics.  A  Wigner-Seitz  cell  method  is  employed  to  identify 
vacancies,  interstitials  and  antisite  defects.  A  lattice  site  with  no  atoms  is  designated  as  a  vacancy, 
the  cell  with  two  atoms  where  one  is  recognized  as  an  interstitial  with  the  site  occupied  by  the 
wrong  atom  type  is  defined  as  an  antisite  defect. 

3.3.  Interatomic  Potentials 


The  interactions  between  Ga-Ga,  Ga-As  and  As-As  are  described  using  analytical  bond-order 
potentials  (BOPs)  developed  by  Albe  et  al.  [23]  The  potentials  have  been  used  to  study  a  wide 
range  of  properties  of  GaAs  compound  structures,  as  well  as  the  pure  phases  of  gallium  and  arsenic, 
including  non-equilibrium  configurations.  The  numerical  results  on  point  defect  properties,  surface 
properties,  and  melting  behavior  have  been  used  to  validate  the  accuracy  and  transferability  of  the 
potentials  [23].  Point-defect  properties  and  surfaces  with  low  As  content  are  found  to  be  in  good 
agreement  with  literature  results,  and  the  molten  GaAs  is  found  to  be  consistent  with  a  structural 
model  based  on  experiments  that  indicate  a  polymerized  arsenic  phase  in  the  melt.  However,  in 
order  to  simulate  the  interaction  of  high  energy  particles  with  GaAs,  the  short-range  interactions 
have  been  further  modified  to  match  the  ‘Universal’  screening  function  Ziegler-Biersack-Littmark 
(ZBL)  potential  [24]  that  gives  a  good  representation  of  high  energy  scattering  of  atoms  in  solids. 
This  is  achieved  by  coupling  the  original  BOPs  and  the  short-range  function  through  the  use  of  the 
Fenni  function 


F(r)  = 


-bf(r-rf  )  , 


(1) 


where  F(r)  quickly  goes  to  one  as  r  increases.  The  total  potential  is  given  by  connecting  the  BOP 
to  the  ZBL  potential: 


V(r)  =  Vm_  l-F(r)  \+VBOpF(  r ) 


(2) 
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which  yields  a  repulsive  potential  that  is  dominated  by  the  ZBL  for  very  short  distance  and  quickly 
approaches  the  BOP  as  r  increases,  as  desired.  The  parameters  of  bf  and  r/  are  adjusted  to  make 
sure  that  the  potential  and  its  first  derivative  are  smoothly  transitioned.  The  values  of  bf  are 
determined  to  be  14,  14  and  7.5  A'1,  and  those  of  r/  are  0.95,  0.95  and  0.65  A  for  Ga-Ga,  Ga-As 
and  As-As  interactions,  respectively.  We  have  employed  these  potentials  to  study  the  threshold 
energies  of  Ga  and  As  [25],  and  the  results  are  in  good  agreement  with  those  obtained  by  both  first 
principles  calculations  and  experimental  methods.  The  fonnation  energies  of  vacancies  are 
determined  to  be  2.0  and  2.4  eV  for  the  As  and  Ga  sublattices,  respectively.  Several  configurations 
of  self-interstitial  atoms  (SIA)  are  found  to  be  stable  in  the  static  lattice,  and  the  most  stable 
interstitial  is  an  As-As<100>  dumbbell  centered  on  As  sites  with  the  formation  energy  of  5.88  eV. 
The  results  are  in  reasonable  agreement  with  those  found  in  ab  initio  calculations.  Nevertheless, 
these  results  suggest  that  the  model  potentials  are  well  suited  for  studying  displacement  cascades 
in  GaAs. 

3.4.  NIEL  Model 

Several  models  have  been  previously  developed  to  determine  the  NIEL  in  Si  [26]  and  GaAs 
[27].  Conventionally,  the  NIEL  can  be  determined  through  the  so  called  Lindhard  energy  partition 
function,  G(T)  [28] 

NIEL  =  a  f Tmax  T—G(T)dT  m 

where  a  is  the  atomic  density  of  the  target  material,  T  the  kinetic  energy  of  the  recoil,  Tmax  the 
maximum  energy  that  can  be  transferred  to  a  recoil  nucleus  by  an  incident  particle  and  Ed  the 
threshold  displacement  energy  which  is  conventionally  taken  to  be  25  eV.  G(T)  is  the  energy 
partition  function  and  do/dT  the  differential  interaction  cross  section.  At  low  energy,  the  recoil 
energy  of  a  PKA,  Epka,  is  close  to  T,  but  at  high  energy  most  of  the  energy  of  a  PKA  contributes 
to  ionization.  Also,  the  number  of  displacements  can  be  estimated  through  G(T)  by  the  Kinchin- 
Pease  model.  Based  on  this  model,  the  NIEL  can  be  rewritten  as 

NIEL  =  a  f  T""L'  F(T)G(T)— dT  (d, 

J  Ed  dT  1  ’ 

where  F(T)  =  EoNf/OA.  It  is  well  known  that  the  number  of  defects  produced  by  a  PKA,  Nf,  is 
linearly  proportional  to  the  recoil  energy  in  the  Kinchin-Pease  model,  which  means  that  in  this 
model  the  NIEL  becomes  proportional  to  the  number  of  defects  produced  by  irradiation.  However, 
it  has  been  recognized  that  the  NIEL  is  not  necessarily  proportional  to  the  number  of  defects,  due 
to  non-linear  processes  that  take  place  in  semiconductors,  particularly  associated  with  the 
formation  of  multiple  disordered  regions  or  amorphous  pockets.  In  the  present  study,  we  have  used 
a  MD  method  to  simulate  defect  production  where  non-linear  processes  will  be  simultaneously 
taken  into  account.  Based  on  our  MD  results,  we  introduce  an  “effective  NIEL”  that  depends  on 
the  number  of  defects  in  the  material,  i.e., 
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(5) 


NIEL  =  a  f  -5 ‘-NmdG(  T)—dT 

JA  0.4  F  dT 

NrMD(T)  is  the  number  of  defects  produced  at  a  given  recoil  energy  of  T,  as  obtained  by  MD 
simulations. 


Many  different  forms  of  screened  Coulomb  potentials  are  proposed  to  calculate  the 
differential  interaction  cross  section,  including  Bohr,  Thomas-Fermi,  Lindhard,  Molere  and  ZBL 
potentials  [29].  Considering  that  the  ZBL  Universal  Formalism  is  based  on  experimentally 
determined  interatomic  potentials,  this  formalism  is  generally  used  to  calculate  the  differential 
interaction  cross  section,  as  given  by 


dcr  _  -na2u  f(  tI/2 ) 
dt  2  tm 


(6) 


where  t  =  t? T/Tmax  is  a  dimensionless  collision  parameter  related  to  recoil  energy  T,  s  is  the  ZBL 
reduced  energy  and  au  is  the  ZBL  Universal  screening  length.  The  function  of f(t1/2)  is  referred  to 
as  the  Thomas-Fermi  scattering  function,  as  detailed  in  references  [26,29].  A  simplified  version  of 
the  Lindhard  energy  partition  function  is  used  for  the  present  study  [26]: 


G(T)= - 

}+Ks(£d) 


(7) 


where  kd,  Sd  and  g( &i)  are  dimensionless  parameters,  and  their  detailed  formulas  and  values  are 
described  elsewhere  [26]. 


4.  RESULTS  AND  DISCUSSIONS 


4.1.  Displacement  Energies 

Within  the  GaAs  microstructure,  a  distinct  threshold  energy  exists  for  each  crystallographic 
direction.  In  the  following  sections,  the  threshold  displacement  energies  on  both  sides  and  for  the 
inner  part  of  unit  stereographic  triangles  will  be  established.  For  a  detailed  description  of  defect 
configurations  and  separation  distance,  we  will  present  these  at  the  end  of  the  recoil  events,  which 
can  be  used  to  explain  the  corresponding  values  of  Ed  in  low  index  directions. 
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Fig.  2  (a)  The  displacement  threshold  energy  around  the  sides  of  the  unit  stereographic  triangles 
for  Ga;  (b)  the  displacement  threshold  energy  around  the  sides  of  the  unit  stereographic  triangles 
for  As. 


4.1.1.  Unit  Stereographic  Triangles 

The  calculated  threshold  energies  for  Ga  and  As  recoil  atoms  on  the  three  borders  of  the  two 
unit  stereographic  triangles  are  illustrated  in  Fig.  2(a)  and  Fig.  2(b),  respectively.  It  was  noticed 
that  both  Ga  and  As  atoms  require  at  least  8  eV  kinetic  energy  to  knock  off  a  neighbor  atom  along 
a  given  direction.  This  is  the  lowest  energy  value  for  which  a  stable  Frenkel  pair  was  created 
through  the  simulations.  For  both  triangles  considered,  the  threshold  energy  of  a  Ga  PKA  varies 
over  a  range  of  8  eV  to  18  eV  along  the  sides  of  the  triangles,  as  shown  in  Fig.  2(a),  while  As  was 
widely  distributed  between  8  eV  and  28  eV  in  Fig.  2(b).  A  lower  value  of  energy  of  8  eV  in  the 
[Ill]  direction  than  that  of  12  eV  in  the  [111]  direction  was  noticed  for  Ga.  The  higher  value  in 
[111]  for  Ga  can  be  explained  by  the  strong  repulsive  energy  between  Ga  and  its  first  nearest 
neighbor  As,  thus  forcing  it  to  return  to  its  equilibrium  site  during  collision.  However,  for  the  As 
case  (Fig.  2b),  the  displacement  threshold  energy  for  these  two  directions  was  8  eV  and  10  eV, 
respectively.  This  could  be  explained  by  the  energy  transfer  efficiency  of  the  heavy  PKA  to  a 
lighter  atom,  thus  the  difference  in  the  two  directions  for  an  As  atom  was  not  significant,  though  a 
larger  distance  to  Ga  in  [111]  direction  would  make  for  a  slight  difference  in  the  value. 
Furthermore,  a  Ga  atom  is  1.13  lighter  than  an  As  atom,  and  thus  it  should  be  easier  for  an  As  to 
knock  off  a  Ga  from  its  lattice  position  than  vice  versa.  However,  in  monatomic  diamond 
structured  semiconductors  (e.g.  Si),  it  has  been  shown  that  Ed  is  equivalent  in  all  (111)  directions 
[30].  The  differences  noticed  here  will  be  further  discussed  in  Section  4.2. 

4.1.2.  Threshold  Displacement  Energy  Contour  Maps 

In  addition  to  the  Ed  determined  along  the  main  crystallographic  directions,  a  large  number 
of  directions  within  the  stereographic  projection  triangles  were  also  taken  into  consideration  as  a 
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function  of  angles.  The  results  are  interpreted  by  Mercator  maps.  As  mentioned  before,  over  three 
thousand  recoil  events  were  perfonned  for  both  Ga  and  As  to  get  a  smooth  and  comprehensive 
distribution,  and  an  in-depth  understanding  of  Ed  dependence  on  crystallographic  orientation,  and 
the  results  are  shown  in  Figs.  3(a)  and  (b)  for  Ga  and  As  recoils,  respectively,  where  the  different 
colors  indicate  different  energy  scales.  Thus,  we  can  consider  the  Ed  as  defining  three  distinctive 
regimes,  namely  zone  1,  2  and  3  as  indicated  in  the  figures.  An  interesting  observation  is  that  there 
is  an  even  higher  value  of  threshold  energy  of  22  eV  for  Ga,  which  occurs  to  trigger  the  recoil 
events  in  the  central  part  of  the  stereographic  projection  triangles,  while  a  maximum  of  18  eV  is 
found  along  the  directions  on  the  sides  for  Ga,  as  described  in  Section  3.1.  This  may  be  due  to  the 
collision  of  the  PKA  with  the  secondary  nearest  neighbor  atoms  after  moving  a  longer  distance. 
Comparing  the  two  triangles  depicting  the  displacement  energies  of  Ga  recoils  (Fig.  3a),  there  is  a 
clear  quantitative  difference  between  the  negative  and  positive  directions.  This  is  especially  true 
in  the  inner  region  (zone  2)  where  the  fraction  of  Ga  recoils  that  produce  persistent  defects  appears 
to  be  higher  in  the  positive  direction,  with  an  initial  energy  of  1 8  eV,  consistent  with  the  larger 
distribution  of  18  eV  on  the  positive  side  than  that  on  the  negative  side.  Threshold  energies  for  Ga 
PKAs  in  the  [111]  direction  and  20  degrees  north  of  the  (110)  pole,  denoted  as  zone  1,  are  4  eV 
greater  than  those  along  the  [111]  direction.  Further,  the  region  enclosed  by  zone  2  reaches  a 
maximum  threshold  of  22  eV,  while  the  corresponding  region  on  the  negative  side  comes  with  a 
maximum  of  18  eV. 
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Fig.  3  (a)  Threshold  energy  map  for  Ga;  (b)  threshold  energy  map  for  As. 

Similarly,  for  As  PKAs  (Fig.  3b),  a  maximum  value  of  28  eV  (Ed)  has  been  found  in  zone  2 
which  can  also  be  interpreted  by  the  same  mechanism  of  collision  with  a  secondary  nearest 
neighbor  as  Ga.  Overall,  while  there  are  small  differences  in  the  region  adjacent  to  the  [100] 
direction  (zone  3),  large  differences  are  noticed  in  zone  1  and  zone  2.  This  indicates  that  there  is  a 
higher  probability  of  producing  a  stable  defect  with  small  energies  in  zone  3,  than  in  either  zone  1 
or  zone  2. 
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Table  1.  Overview  of  displacement  energies  compared  with  experiment  data. 

(Ed.ave,  average  displacement  energy  over  three  thousand  directions  for  each  type;  Ndirection, 
_ _ _ number  of  directions  simulated  for  each  type.) _ _ 


Recoil 

Ed(eV) 

type 

Ndirection 

Ed,ave 

(eV) 

[100] 

[111] 

[110] 

[  1 10 

] 

[111] 

[111] 

Experiment 

Ga 

1828 

13.38 

10 

10 

14 

14 

12 

8 

10.0+0.7  [32] 

As 

1828 

12.75 

10 

10 

18 

18 

10 

8 

9  [32] 

4.1.3.  Overview  of  Ed 

The  average  displacement  threshold  energies  in  two  oppositely  oriented  stereographic  triangles  for 
each  atom  type  are  summarized  in  Table  1.  The  average  threshold  energies  over  the  statistically 
significant  sample  size  (few  thousand  directions)  combining  both  positive  and  negative  directions 
were  found  to  be  12.75  and  13.38  eV  for  As  and  Ga,  respectively.  The  almost  identical  value  of 
Ed.ave  (small  difference  of  0.63  eV)  is  to  be  expected  due  to  the  many  physical  similarities  between 
these  two  atoms.  Moreover,  there  have  also  been  several  experimental  studies  in  GaAs  on 
examining  primary  damage  response.  The  lowest  energy  threshold  from  the  experimental 
estimation  of  10.0+0.7  eV  and  9  eV  [31]  is  consistent  with  8  eV  at  the  minimum  for  both  species 
from  our  simulations.  This  indicates  a  good  match  with  the  established  data  in  crystallographic 
non-equivalent  crystal  directions.  In  addition,  since  the  effective  displacement  energy  is  typically 
about  50%  higher  than  the  threshold  [32]  an  estimation  value  of  15  eV  on  top  of  experiment  data 
is  consistent  with  what  we  have  obtained  through  MD  simulations  (13. 1+0. 3eV).  Hence,  they  are 
in  satisfactory  agreement. 

4.2.  Displacement  Cascades  and  Effective  NIEL 
4.2.1.  Defect  Generation 

Similar  to  other  semiconductors,  the  displacement  cascades  modeled  in  the  present  work 
exhibit  two  basic  phases,  i.e.  a  ballistic  phase  and  a  recombination  phase.  Here,  the  former  lasts  a 
few  tenths  of  ps,  during  which  the  primary  recoil  atom  initiates  an  avalanche  of  atomic 
displacements,  resulting  in  many  atoms  being  displaced.  After  the  number  of  displacements 
reaches  a  maximum,  the  system  starts  to  cool  down  and  the  follow-up  recombination  occurs, 
causing  some  of  these  displacements  to  recombine  with  the  empty  lattice  sites.  At  about  5  ps,  the 
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cascade  process  reaches  a  stable  state  with  occasional  interstitial  rotation,  but  without  significant 
defect  migration.  The  peak  and  final  damage  states  of  a  typical  10  keV  cascade  in  GaAs  are  shown 
in  Figs.  4(a)  and  (b),  respectively,  where  only  displaced  atoms  are  plotted.  The  interstitials  are 
represented  by  large  spheres  and  vacancies  by  small  ones,  as  indicated  in  the  legend.  The  Ga  PKA 
in  this  simulation  traveled  about  10  nm  before  it  came  to  rest  in  the  lattice,  generating  multiple 
damage  domains  along  its  path.  These  domains  form  distinct  regions  separated  from  each  other, 
within  which  a  large  number  of  atoms  are  displaced.  The  final  damage  is  very  similar  to  the  peak 
damage  in  the  case  of  GaAs,  but  with  slightly  smaller  number  of  displacements.  It  is  surprising 
that  the  recombination  of  displaced  atoms  with  the  empty  lattice  sites  is  small,  in  contrast  to  the 
displacement  cascades  in  SiC  [33]  and  metals  [34]  where  50%  and  90%  of  the  displaced  atoms 
recombine  with  vacancies,  respectively.  This  may  be  associated  with  the  low  threshold  energies 
of  GaAs  [25]  where  a  PKA  is  able  to  create  large  disordered  regions  along  its  path,  leading  to  the 
formation  of  amorphous  domains,  which  prevents  defect  recombination  at  low  temperatures.  Fig. 
4(b)  clearly  shows  that  the  damage  consists  of  multiple  disordered  regions  or  defect  clusters,  in 
contrast  to  SiC  where  most  defects  are  single  interstitials  and  mono-vacancies,  and  only  a  small 
proportion  of  the  interstitials  are  found  in  clusters. 


•  Gaj  •  Ga^ 


Fig.  4  Atomistic  configurations  of  a  typical  10  keV  cascade  in  GaAs  at  (a)  the  peak  damage  and 
(b)  the  final  damage  state,  where  only  defects  are  shown.  The  defects  are  identified  by  size  and 
color,  as  indicated  by  different  spheres. 
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To  better  understand  the  defect  production  in  GaAs,  the  numbers  of  antisite  defects  and 
atoms  displaced  into  interstitial  positions  are  shown  in  Fig.  5  as  a  function  of  time  for  the  10  keV 
cascade  demonstrated  above,  where  Ga  and  As  components  are  plotted  separately.  The  increase  in 
the  number  of  displaced  atoms,  Nd,  during  the  collisional  stage  reaches  a  maximum  at  tpeak~  0.5 
ps,  after  which  Nd  decreases  in  a  recombination  phase.  As  noted  before,  the  number  of  displaced 
defects  decreases  only  slightly  to  final  values  within  1.3  ps.  However,  the  recombination  of 
displaced  atoms  with  vacancies  is  small  within  the  relaxation  phase,  as  shown  in  Figs.  4(a)  and 
(b),  and  only  18%  of  the  displaced  atoms  recombine  with  vacancies.  Similar  to  other 
semiconductor  compounds,  most  antisite  defects  are  generated  during  the  collisional  phase,  when 
the  maximum  in  Nd  is  reached,  and  the  increase  in  antisite  defects  is  very  small  during  relaxation. 
This  may  suggest  that  the  formation  of  antisite  defects  is  due  to  replacement  collision  sequences, 
although  they  are  very  short.  However,  the  cascade  lifetime,  defined  as  the  time  during  which  the 
variation  in  the  final  number  of  defects  can  be  neglected,  is  very  short  as  ~  1.2  ps,  and  the  defects 
are  quickly  quenched,  thus  forming  high  disordered  domains.  Previously,  several  models  have 
been  proposed  to  explore  the  mechanisms  associated  with  amorphization  due  to  ion  irradiation 
[34],  including  direct-impact  amorphization  within  an  individual  collision  cascade,  the 
accumulation  of  defects  due  to  cascade  overlap,  defect-stimulated  nucleation,  and  growth 
processes,  or  a  combination  of  these  processes.  The  observed  generation  of  amorphous  pockets  in 
GaAs,  which  are  formed  locally  within  a  cascade,  is  likely  to  be  associated  with  direct-impact 
amorphization  in  the  model.  Another  important  feature  is  that  the  number  of  Ga  defects  is  roughly 
the  same  as  that  for  As,  which  is  expected  because  of  the  similar  threshold  energies  for  Ga  and  As 
in  GaAs,  i.e  the  threshold  displacement  energies  are  13.75  eV  for  Ga  and  14.3  eV  for  As  averaged 
over  more  than  3600  simulated  recoil  events  [25]. 


Time  (ps) 


Fig.  5  The  number  of  displaced  atoms,  including  Ga  and  As  displacements  and  antisites,  as  a 
function  of  time  for  the  10  keV  cascade  shown  in  Fig.  4. 
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Figure  6  shows  the  final  defect  distributions  and  cascade  morphologies  produced  by  2,  5  and 
20  keV  cascades  in  GaAs.  The  defects  generated  by  a  PKA  with  an  energy  less  than  2  keV 
distribute  on  one  cascade  region,  and  consist  of  interstitials  and  vacancies,  rather  than  antisite 
defects,  as  demonstrated  by  Fig.  6(a).  However,  a  PKA  generates  multiple  sub-branches  of  the 
cascade  along  its  path  when  its  energy  is  larger  than  2  keV,  which  is  clearly  displayed  in  Figs.  6(b) 
and  (c),  forming  distinct  regions  separated  from  each  other  with  a  high  energy  density  in  the  cores 
of  the  sub-branches.  For  example,  a  5  keV  PKA  in  Fig.  6(b)  creates  two  distinct  regions  each  of 
which  has  its  own  core  of  high  defect  density.  This  suggests  that  the  nonlinear  effects  take  place 
in  high  energy  cascades  in  GaAs,  which  may  be  important  in  calculating  the  NIEL  with  effects  not 
captured  in  the  simplified  Kinchin-Pease  model. 


2  keV 


Fig.  6  The  final  defect  states  of  (a)  2  keV,  (b)  5  keV  and  (c)  20  keV  cascade  in  GaAs  from  which 
can  be  clearly  seen  the  transition  from  a  single  pocket  of  atomic  displacements  to  multiple 
subcascades. 


The  final  number  of  defects  created  in  displacement  cascades  is  an  important  outcome  of 
MD  simulations  because  of  its  applications  to  theories  and  models  of  radiation  damage  evolution. 
This  parameter  can  be  easily  determined  from  MD  simulations  and  compared  with  the  standard 
formula  of  Norgett,  Robison  and  Torrens  [35]  defined  by 


N  =0.8E  /  2E 

NRT  PKA  d 


(8) 


where  Epka  is  the  PKA  energy  that  is  elastically  deposited  into  the  target,  and  Ed  is  the  mean 
threshold  energy  that  is  taken  to  be  14  eV  [25].  Data  for  the  final  number  of  Frenkel  pairs,  Nf,  at 
the  end  of  the  displacement  cascades  are  shown  as  a  function  of  PKA  energy  in  Fig.  7,  along  with 
the  NRT  values  for  comparison.  The  number  of  Ga  interstitials,  As  interstitials  and  antisite  defects 
are  presented  separately,  and  the  points  denote  the  mean  values  for  the  displacement  cascades 
investigated  at  the  given  energy.  It  is  of  interest  to  note  that  the  total  number  of  defects  produced 
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in  GaAs  is  generally  much  higher  than  that  predicted  by  the  NRT  formula,  in  contrast  to  SiC  where 
the  total  number  of  the  defects  is  smaller  than  the  NRT  value  [36].  This  may  be  associated  with 
the  direct-impact  amorphization  that  occurs  in  GaAs,  which  prevents  the  fast  annihilation  of 
interstitials  with  vacancies  during  the  relaxation  phase.  The  defect  production  efficiency,  defined 
as  the  ratio  of  Nf  to  Nnrt,  ranges  from  3.0  to  5.0  and  increases  with  increasing  PKA  energy,  rather 
than  decreasing  with  increasing  PKA  energy  as  revealed  in  SiC  [36]  and  metals  [34].  This  suggests 
that  defect  production  efficiency  in  GaAs  is  higher  than  that  found  in  metals  or  alloys  because  of 
the  absence  of  a  thermal  spike  in  GaAs.  It  is  also  clear  that  there  exist  two  distinct  regions  of  defect 
production,  one  exhibiting  the  linear  behavior  of  cascades  when  PKA  energy  is  less  than  2  keV 
and  a  non-linear  region  with  enhanced  defect  production  related  to  the  formation  of  amorphous 
domains  for  PKA  energies  larger  than  2  keV,  as  shown  in  the  insert  of  Fig.  7.  This  crossover  in 
defect  production  may  account  for  the  different  mechanisms  for  energy  partition  in  GaAs.  The 
NRT  values  are  imposed  for  comparison.  Data  points  are  the  mean  values  at  each  energy. 


Fig.  7  The  number  of  surviving  defects  produced  in  the  cascades  as  a  function  of  PKA  energy, 
where  Ga  and  As  interstitials,  as  well  as  antisite  defects,  are  presented  separately. 
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EPKA/Npk  (eV/atom) 


Fig.  8  The  total  defects  of  the  final  damages  as  a  function  of  the  energy  density  deposited  by  a 
PKA,  where  the  solid  line  represents  the  fits  to  MD  data. 


4.2.2.  Effective  NIEL 

Based  on  the  MD  simulations,  the  number  of  Frenkel  pairs  produced  can  be  extrapolated  into 
higher  energies,  and  thus  it  is  possible  to  develop  a  relationship  describing  defect  production  and 
PKA  energy.  Previous  simulations  of  displacement  cascades  in  a  number  of  materials  have 
revealed  that  the  number  of  Frenkel  pairs  displays  a  power-law  dependence  on  PKA  energy, 
weakly  depending  on  materials  and  temperature.  If  we  assume  that  the  energy  transfer  to  atoms 
occurs  just  before  the  damage  peak  (maximum),  the  energy  density  deposited  by  a  PKA  is  equal 
to  ped=EpKA/Npk,  where  Epka  is  the  PKA  energy  and  NPk  is  the  number  of  displacements  generated 
at  the  peak  damage.  The  number  of  surviving  defects  (Frankel  pairs  produced)  can  be  described 
by 


ME„Jp' 


ed 


(9) 


where  the  PKA  energy,  Epka,  is  given  in  keV.  The  dimensionless  parameters  A,  n  and  m  can  be 
determined  by  fitting  to  the  MD  data,  and  the  values  for  GaAs  are  approximately  7635.95,  2. 1499 
and  1 . 1 108,  respectively.  The  total  number  of  defects,  that  includes  interstitials  and  antisites, 

is  plotted  as  a  function  of  energy  density  in  Fig.  8,  where  the  solid  line  indicates  the  values 
calculated  by  Eq.  (9).  It  can  be  seen  that  our  analytical  model  is  relatively  close  to  the  MD 
simulations.  NrMI)  decreases  with  increasing  energy  density,  which  is  expected  because  a  large 
number  of  atoms  are  displaced  from  their  sites  at  higher  PKA  energies,  giving  rise  to  the  lower 
energy  density.  Practically,  a  rate  of  degradation  in  electronic  devices  can  be  described  by  the 
number  of  defects  in  the  materials,  but  experimentally,  degradations  are  represented  by  the  NIEL 
which  is  assumed  to  be  proportional  to  the  number  of  defects  produced  by  irradiation.  Combining 
Eqs.  (5)  and  (9),  the  effective  NIEL  in  GaAs  is  determined  for  incident  proton,  alpha  and  Xe 
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particles,  and  the  results  are  shown  in  Fig.  9,  along  with  other  calculations  for  comparison.  The 
striking  fact  is  that  the  general  behavior  is  very  similar  in  all  of  these  calculations,  i.e.  the  NIEL 
increases  with  increasing  incident  energy  until  the  peak  NIEL  and  then  decreases.  Also,  the  peak 
position  and  magnitude  of  the  NIEL  strongly  depend  on  the  incident  particles,  and  the  peak  NIEL 
shifts  to  higher  energies  with  increasing  mass  of  incident  particles.  For  instance,  the  peak  NIEL 
for  Xe  is  more  than  5000  times  higher  than  that  for  protons,  and  the  peak  occurs  at  an  energy  that 
is  about  1000  times  higher  than  that  for  protons.  Messenger  et  al.  used  a  similar  model  to  calculate 
the  NIEL  in  numerous  materials  [26],  including  GaAs,  but  assumed  that  displacement  damage  is 
proportional  to  the  kinetic  energy  of  the  recoil.  Their  results  are  also  included  in  Fig.  9,  as  indicated 
by  the  dotted  lines.  Although  the  general  tendency  of  the  NIEL  is  similar,  there  appears  to  be 
distinguished  differences  at  both  low  and  high  energies,  which  may  point  to  fundamentally 
different  mechanisms  of  displacement  damage  in  GaAs.  As  described  before,  a  PKA  generates  a 
single  pocket  of  atomic  displacements  when  its  energy  is  below  2  keV,  or  single  point  defects 
when  its  energy  is  less  than  500  eV,  where  the  damage  is  in  the  linear  regime  of  defect  generation 
and  is  proportional  to  the  recoil  energy.  This  may  provide  an  explanation  for  why  the  NIEL  at  low 
energies  is  in  reasonable  agreement  with  that  calculated  by  Messenger  et  al  [26].  For  higher 
energies,  the  sub-branches  or  subcascades  generated  by  recoil  nuclei  must  be  taken  into  account, 
which  describes  the  non-linear  behavior  of  defect  production  in  GaAs  and  cannot  be  predicted  by 
the  simple  Kinchin-Pease  model.  For  instance,  the  number  of  atomic  displacements  by  a  5  keV 
PKA  in  GaAs  is,  on  average,  around  421,  that  is  about  3  times  the  number  of  displacements 
predicted  by  the  Kinchin-Pease  model.  The  multiple  amorphous  pockets  fonned  within  a  cascade 
give  rise  to  the  number  of  defects  produced,  resulting  in  larger  NIEL  at  high  recoil  energies  than 
that  determined  by  Messenger  et  al  [26].  In  addition,  the  NIEL  determined  using  the 
simulation  code  based  on  transport  of  ions  in  matter  transport  of  ions  in  matter  (TRIM)  code 
for  incident  protons  from  0.2  to  2  MeV  in  GaAs  [14]  is  added  in  Fig.  9  for  comparison.  The 
data  points  are  very  close  to  the  NIEL  of  Messenger  et  al.  and  smaller  than  the  calculations 
based  on  MD  simulations,  which  is  expected.  In  the  TRIM  calculation,  the  vacancy  production 
rate  as  a  function  of  position  was  combined  with  the  total  energy  loss  data  to  determine  the 
vacancy  production  rate  as  a  function  of  proton  energy.  It  should  be  noted  that  TRIM 
calculations  do  not  take  defect  kinetics  (e.g.  recombination)  into  account.  The  sub-cascade 
mechanisms  outlined  herein  seem  to  overcompensate  for  any  point  defect  kinetics  ignored 
by  TRIM.  Since  the  vacancy  concentration  was  converted  to  damage  energy  using  the 
modified  Kinchin-Pease  approximation  to  yield  the  NIEL  as  a  function  of  proton  energy,  it  is  not 
surprising  that  the  NIELs  calculated  by  Messenger  et  al.  and  TRIM  are  in  good  agreement. 
Nevertheless,  we  have  developed  an  MD-based  NIEL  model  to  determine  the  NIEL  in  GaAs, 
which  can  be  applied  to  a  variety  of  incident  particles.  It  is  also  suggested  that  further 
experimental  measurements  are  needed  for  higher  fidelity  comparisons  with  this  and 
other  methods. 


15 


Approved  for  public  release;  distribution  is  unlimited. 


lO"  10' 

Energy  (MeV) 


Fig.  9  The  NIEL  obtained  based  on  the  MD  results  for  several  different  incident  particles  in  GaAs, 
where  the  values  of  the  NIEL  calculated  by  an  analytical  approach  [26]  are  imposed  for 
comparison. 


5.  CONCLUSIONS 

The  threshold  displacement  energies  for  Ga  and  As  atoms  were  investigated  using  a 
molecular  dynamics  method  with  a  bond-order  potential  for  GaAs,  and  the  results  demonstrated  a 
weak  dependence  on  the  crystallographic  direction  in  which  the  atom  is  displaced.  The  different 
types  of  lattice-atom  vibrational  modes  and  atomic  sizes  are  believed  to  account  for  this  slight 
dependence.  By  considering  the  crystallographic  orientations  of  GaAs,  a  total  of  3600  recoils  are 
simulated  to  study  the  special  distribution  of  threshold  energies  in  four  stereographic  triangles. 
The  threshold  energy  for  As  recoils  ranged  from  8  to  28  eV,  while  it  ranged  from  8  to  22  eV  for 
Ga  recoils.  The  average  threshold  along  the  directions  whose  crystallographic  projections  are 
located  on  the  central  part  of  a  triangle  is  generally  higher  than  on  its  sides.  Nonetheless,  it  has 
also  been  shown  that  the  average  displacement  energies  were  similar  for  both  atoms,  which  is 
reasonable  since  the  two  atoms  have  similar  physical  properties  like  mass  and  short-range 
repulsion.  Moreover,  for  low-energy  displacements,  the  defect  types  generated  by  recoil  events 
along  low-index  crystallographic  directions  are  mostly  in  a  simple  form  of  Frenkel  pairs.  The 
separation  distances  of  Frenkel  pairs  in  these  directions  are  mostly  around  4.8  A  at  an  average  of 
4.69  A  except  those  along  the  [100]  direction,  which  become  shorter  with  a  small  threshold  energy 
of  10  eV  and  a  closer  distance  from  the  initial  position  to  a  stable  octahedral  site.  However, 
complex  defect  configurations  along  some  high-index  directions  are  observed,  which  results  in 
higher  threshold  energies. 
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Also,  we  have  employed  large-scale  molecular  dynamics  (MD)  simulations  to  study  defect 
production  and  its  detailed  evolution  for  PKA  energies  from  500  eV  to  20  keV  in  GaAs.  At  low 
energies,  most  surviving  defects  are  single  interstitials  and  vacancies,  and  only  20%  of  the 
interstitial  population  is  contained  in  clusters.  A  PKA  with  an  energy  less  than  2  keV  generates  a 
single  pocket  of  atomic  displacements,  which  is  associated  with  the  linear  regime  of  displacement 
damage.  Direct-impact  amorphization  occurs  with  a  high  degree  of  probability  when  PKA  energies 
are  larger  than  2  keV,  which  involves  non-linear  defect  production.  The  transition  from  linear  to 
non-linear  defect  production  plays  an  important  role  in  the  calculation  of  the  NIEL  in  GaAs.  Based 
on  the  MD  results,  we  have  further  detennined  the  damage  density  within  a  cascade  core  that  is 
used  to  develop  a  model  for  calculating  the  NIEL  in  GaAs,  which  can  then  be  used  to  predict  the 
displacement  damage  degradation  induced  by  space  radiation  on  electronic  components. 

The  calculated  NIEL  results  are  compared  with  the  available  data  detennined  by  different 
models.  The  agreement  between  our  calculations  and  other  values  of  the  NIEL  for  protons  is 
reasonably  good  at  low  energies,  but  there  is  a  discrepancy  at  higher  energies,  which  is  expected 
because  of  non-linear  effects  neglected  by  the  Kinchin-Pease  model.  The  model  developed  in  the 
present  study  provides  a  general  method  for  accurately  calculating  the  NIEL  for  fast  neutrons, 
protons  and  other  heavy  particles,  and  can  also  be  applied  to  other  semiconductor  materials. 
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LIST  OF  ACRONYMS 


As 

Arsenide 

BOPs 

Bond-order  potentials 

eV 

Electron  volt 

GaAs 

Gallium  arsenide 

InP 

Indium  phoshide 

keV 

1000  eV 

MD 

Molecular  dynamics 

MOLDY 

MOLecular  Dynamics  simulation  code 

NF 

the  final  number  of  Frenkel  pairs 

NIEL 

Non-ionizing  energy  loss 

PKA 

Primary  knock-on  atom 

Ps 

Pico- second 

Rf 

Fitting  parameter 

SIA 

Self-interstitial  atoms 

SiC 

Silicon  carbide 

Si 

Silicon 

SREC 

Space  radiation  effects  community 

TRIM 

Transport  of  ions  in  matter 

ZBL 

Ziegler-Biersack-Littmark 

21 


Approved  for  public  release;  distribution  is  unlimited. 


DISTRIBUTION  LIST 


DTIC/OCP 

8725  John  J.  Kingman  Rd,  Suite  0944 
Ft  Belvoir,  VA  22060-6218  1  cy 

AFRL/RVIL 

Kirtland  AFB,  NM  87117-5776  2  cys 

Official  Record  Copy 

AFRL/RVSW/David  Cardimona  1  cy 


22 


Approved  for  public  release;  distribution  is  unlimited. 


